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Abstract — In relay-aided wireless transmission systems, one of 
the key issues is how to manage the energy resource at the source 
and each individual relay, to optimize a certain performance 
metric. This paper addresses the sum rate maximized resource 
allocation (RA) problem in an orthogonal frequency division 
modulation (OFDM) transmission system assisted by multiple 
decode-and-forward (DF) relays, subject to the individual sum 
power constraints of the source and the relays. In particular, 
the transmission at each subcarrier can be in either the direct 
mode without any relay assisting, or the relay-aided mode with 
one or several relays assisting. We propose two RA algorithms 
which optimize the assignment of transmission mode and source 
power for every subcarrier, as well as the assisting relays and the 
power allocation to them for every relay-aided subcarrier. First, 
it is shown that the considered RA problem has zero Lagrangian 
duality gap when there is a big number of subcarriers. In this 
case, a duality based algorithm that finds a globally optimum RA 
is developed. Most interestingly, the sensitivity analysis in convex 
optimization theory is used to derive a closed-form optimum 
solution to a related convex optimization problem, for which the 
method based on the Karush-Kuhn-Tucker (KKT) conditions 
is not applicable. Second, a coordinate-ascent based iterative 
algorithm, which finds a suboptimum RA but is always applicable 
regardless of the duality gap of the RA problem, is developed. The 
effectiveness of these algorithms has been illustrated by numerical 
experiments. 

Index Terms — Orthogonal frequency division modulation, re- 
source allocation, decode and forward, relaying, Lagrangian 
duality gap, dual decomposition method, energy efficiency. 



I. Introduction 

Relay-aided cooperative transmission finds plenty of 
promising applications when it is difficult to install multiple 
antennas at the same radio equipment, and therefore has been 
attracting intensive research interest in both academia and in- 
dustry lately JTJ. Low complexity yet efficient protocols, such 
as amplify and forward (AF) as well as decode and forward 
(DF), have been proposed to simplify the implementation with 
practical devices |2), (3). Typically, both protocols propose to 
carry out a relay-aided transmission within two time slots, 
namely a broadcasting slot and a relaying slot. In J3J, the 
AF/DF which fixes every transmission in the relay-aided mode 
independently of source-relay channel conditions is referred 
to as fixed relaying AF/DF. In fact, improved protocols may 
be built for better performance. For instance, selection re- 
laying AF/DF, which selects either the direct or relay-aided 
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transmission mode depending on channel conditions, has been 
proposed to improve spectral efficiency JJJ. In particular, the 
direct mode, which refers to the direct source to destination 
transmission without any relay assisting, is used when the 
source-destination channel gain is higher than the source- 
relay channel gain. Most interestingly, not fixed but selection 
relaying DF can achieve full diversity (3J. 

We consider in this paper a point to point orthogonal 
frequency division modulation (OFDM) transmission system 
aided by multiple DF relays. The motivation behind this is 
that the OFDM transmission has been widely recognized for 
current and future wireless systems, thanks to its flexibility to 
incorporate dynamic resource allocation (RA) for performance 
improvement (HJ. In such a transmission system, one of the key 
issues is how to decide for every subcarrier the transmission 
mode and assisting relays if the relay-aided mode is chosen, 
and the power of the source and every individual relay, to 
maximize a certain objective related to system performance. 
Obviously, this RA problem is more complicated compared 
to those for conventional OFDM systems without relays. 
Therefore, RA algorithms are solicited for relay-aided OFDM 
systems. 

To date, some related research works have been reported. 
To name a few, RA algorithms have been proposed in ||5l- 
(8) for OFDM systems aided by AF relays, and in ll9l-lfT4l 
for multi-user OFDM systems aided by DF relays. As for 
OFDM systems aided by DF relays, RA algorithms have 
been proposed in lfT31 to minimize sum power under rate 
constraints, and in |[T6l . IfTTl to maximize sum rate subject 
to power constraints, when a single relay exists. However, 
at a subcarrier in the direct mode the source is idle during 
the relaying slot, which wastes spectrum resource. To address 
this issue, rate-optimized RA algorithms, which allow for the 
source to destination transmission during the two slots at the 
subcarrier in the direct mode, have been proposed in [HB3-GT]]. 

So far, the majority of proposed RA algorithms, as the 
aforementioned ones, restrict that at most one relay can assist 
the source at every relay-aided subcarrier. In fact, when there 
are multiple relays available, allowing not just one, but each 
of them to be eligible for assisting can exploit all degrees of 
freedom in the system to improve performance. For illustration 
purposes, example patterns of selecting single or multiple 
assisting relays are shown in Figure Q] For example, it has 
been shown in l22l that the sum power can be reduced in 
a multi-user OFDM system if multiple DF relays assist the 
transmission at every relay-aided subcarrier. 

In this paper, we address the sum rate maximized RA prob- 
lem in an OFDM system aided by multiple DF relays subject 
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Fig. 1. Example patterns of selecting (a) one relay or (b) multiple relays 
assisting at each relay-aided subcarrier, where a black block at the i-th row 
and j-th column indicates the j'-th relay rj assists the transmission at the i-th 
subcarrier Sj. 

to the individual (per device, i.e. a source or a relay) sum 
(over all subcarriers) power constraints^ of the source and the 
relays. In particular, one or several relays may cooperate with 
the source to transmit at every relay-aided subcarrier. When the 
sum power consumed by the source and that by every relay are 
fixed, the optimum RA to this problem leads to the maximum 
energy efficiency for the system under consideration, because 
the total energy consumed by the source and the relays for 
transmitting per information bit is minimized. Specifically, our 
contributions lie in the following aspects: 

> when a big number of subcarriers is used, it is shown that 
the duality gap of the RA problem is equal to zero, based 
on the same idea first proposed in 11231 . Assuming the 
number of subcarriers is sufficiently large, we develop a 
duality based algorithm which finds a globally optimum 
RA for the considered problem. Most interestingly, the 
sensitivity analysis in convex optimization theory is used 
to derive a closed-form optimum solution to a related 
convex optimization problem, for which the method based 
on the Karush-Kuhn-Tucker (KKT) conditions is not 
applicable. 

« we develop a coordinate-ascent based iterative algorithm, 
which finds a suboptimum RA but is always applicable 
regardless of the duality gap of the RA problem. Specif- 
ically, this algorithm produces a successive set of RAs 
with nondecreasing sum rate until convergence. 
The remainder of this paper is organized as follows. In 
the next section, the OFDM system under consideration is 
described and the RA problem is formulated. Then, the duality 
based algorithm and the iterative algorithm are developed in 
Sections III and IV, respectively. In Section V, the effective- 
ness of the proposed algorithms is illustrated by numerical 
experiments. Finally, some conclusions wrap up this paper in 
Section VI. 

II. System description and RA problem 

FORMULATION 

A. System description 

We consider the OFDM transmission from a source to a 
destination aided by N DF relays collected in the set = 

1 In this paper, the word "individual" will mean per device and "sum" will 
refer to a summation over all subcarriers, unless otherwise stated. 
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TABLE I 

Channel coefficient at subcarrier k between any two of the 

SOURCE, Ti, AND THE DESTINATION. 



source to destination 


source to r 4 


r, to destination 









= 1, • • ■ ,N}. All links are assumed to be frequency 
selective, and OFDM with properly designed cyclic prefix is 
used to transform every link into K parallel channels, each at 
a different subcarrier facing flat fading. At every subcarrier, 
the transmission of a symbol is in either the direct mode, or 
the relay-aided mode spanning across two equal-duration time 
slots, namely the broadcasting slot and the relaying slot. We 
assume the destination decodes the signal samples received 
at each subcarrier separately from those received at any other 
subcarrier. 

We make the following assumptions about the RA in the 
system. First, the RA is determined by an algorithm running 
at a central controller, which knows precisely the noise power 
at each node, as well as the channel coefficients at every 
subcarrier from the source to every r^, from the source to the 
destination, and from every to the destination, respectively. 
Second, all channels remain unchanged within a sufficiently 
long duration, over which RA can be carried out accordingly. 
Third, the RA information can be reliably disseminated to the 
source, every relay, and the destination. 

Let's consider the transmission of a unit-variance symbol 8 
at subcarrier k. The coefficient of the channel between any two 
of the source, r,, and the destination, are notated according 
to Table Q] We first describe the transmission in the relay- 
aided mode. The source first emits in the broadcasting slot 
the symbol \J P s P s ,kfi as illustrated in Figure [2] a, where P s 
and P S) f. represent the source sum power and the fraction of 
that sum power allocated to the transmission at subcarrier 
k, respectively. At the end of this slot, both the destination 
and the relays receive the source signal. The signal samples 
received at the destination and can be expressed respectively 
by 

Vk = ^PsPsMhsAk) + n k (1) 

and 

y ri ,k = V PsP S Mh SiIt {k)0 + n rtik , (2) 

where nu and n ri ^ represent the corruption of the additive 
white Gaussian noise (AWGN) at the destination and r^, 
respectively. We assume Vi, n ri _k is a zero-mean circularly 
Gaussian random variable with the variance a 2 . The signal to 
noise ratio (SNR) at r,; can be computed as P s .kG s , ri (k), where 

G s ,r;(fc) = P3 ^ 3 a 'j ^ represents the normalized channel 
power gain from the source to r^. 

For subcarrier k, we define ^ as the set containing all 
relays sorted in the increasing order of G Sjli (fc), r.;(fc) as the 
i-th relay in ^k, and $k(i) as the set containing relays in 
$k with indices from i to N. We assume the r,; that has the 
minimum G s . ri (k) among all assisting relays for subcarrier k 
is ib k (k) where bk represents the index of that r; in ^J/^, and 
all rj's in v &fe(&fe) decode the received samples to recover 9. 



TO APPEAR IN IEEE JOURNAL ON SELECTED AREAS IN COMMUNICATIONS 



3 



After the recovery of 9, all relays in ^ k (b k ) transmit 
simultaneously to the destination in the relaying slot, which 
in effect establishes a distributed multiple input and single 
output (MISO) transmission link as illustrated in Figure |2]b. 
Specifically, r$ g *&k(b k ) transmits w Ti 9, where w Ti is 
the complex weight for transmit beamforming and satisfies 
\w ri \ 2 = P Ti Pri,k- Pii and P Ti ,k represent the sum power of 
Ti and the fraction of that sum power allocated to subcarrier k, 
respectively. To have the relays' signals add coherently when 
received at the destination, w ri = y/P^P^.k^ 1 arg ( - h ''- d( - k ^ 
is used, where arg (h ru d{k)) stands for the phase of /i ri .d(fc)- 
It should be noted that this transmission protocol enables 
a flexible use of all relays opportunistically through a gen- 
eral form of adaptive transmit beamforming, in that b k and 
{P Iu kVi G 'J'fc(frfc)} are determined by the RA algorithm 
depending on channel state information, as will be developed 
later. At the end of the relaying slot, the signal sample received 
at the destination is denoted by 

Zk= V p r, p u.k\KAk)\0 + v k , (3) 

where v k represents the AWGN corruption at the destination. 
We assume rik and v k are independent zero-mean circularly 
Gaussian random variables with the same variance a 2 . 

Finally, y k and z k are processed at the destination. From 
([T]i and (0, it can be seen that the relay-aided transmission 
at subcarrier k is in effect over a channel with a single input 
and two outputs. To achieve the capacity of this channel, the 
maximum ratio combining (MRC) should be used, i.e., the 
destination combines y k and Zf. to construct a decision variable 

Ck =V ' P s P s ,k{hs.d{k))* : y k + 

( V P ^Pr„k\KAk)\)*z k , (4) 

r«e* fc (&jb) 

which is then decoded (please refer to page 179 of ll24l for 
more details). After mathematical arrangements, the SNR for 
decoding Ck is derived as 

Vk = Ps.kG sA (k) +( ]T Pr ilk G rud (k)) (5) 

where G s , d (fc) = ^fiJ! and G n , d (k) = 

represent the normalized channel power gains from the source 

to the destination and from r; to the destination, respectively. 




(a) In the broadcasting slot (b) In the relaying slot 



Fig. 2. Illustration of a relay-aided transmission where T2 and r3 assist 
relaying at subcarrier k. 



To ensure the reliable recovery of 8 at every relay in ^k{b k ), 
the source transmission rate at subcarrier k should not be 
higher than log 2 (l + P s .kG slb (k)(k)) bits/two-slots (bpts). 
Moreover, the source transmission rate should not be higher 
than log 2 (l + r/k) bpts to ensure reliable decoding of Ck 
at the destination. Therefore, the source transmission rate at 
subcarrier k in the relay-aided mode should be J3] 

R k>1 =min(log 2 (l + P s ^G s ^j fc) (fc)),log 2 (l + ?7 fe )) 
= log 2 (l + min(?7 fc , P Stk G StIb j k) (k))) bpts. (6) 

When bk is fixed, R k ,i is a concave function of 
{P s ,k, Pii.k\y T i € due to the following reasons. 

First, r/k and P s .kG S7lbk (k){k) are both concave 
functions of {P s ,k, Prt,k\V r i G an d therefore 

minjr^, P s ,fcG s . ri)fc (fc)(fc)} is a concave function of 
{P s ,k, Pii.k\y T i G because the minimum of two 

concave functions is still a concave function. Second, 
log 2 (l + x) is an increasing and concave function of x. 
Therefore, Rk,i, as the composition of log 2 (l + x) and 
x = min{?7fc, P s ,kG SIb (k)(k)}, is a concave function 
of {P s .k, -Pri.fclVr,; <G when bk is fixed, because the 
composition of an increasing concave function and a concave 
function is still a concave function l25l . 

As for the direct transmission mode at subcarrier k, we 
consider a more efficient protocol compared to restricting the 
source to transmit only in the broadcasting slot as in the 
related works J3], |[T5l - lfT7l . Specifically, the source emits 
two independent symbols in the two slots, respectively, and 
only the destination decodes the corresponding two received 
signal samples. We assume the AWGN corruptions for the 
two received samples are independent zero-mean circularly 
Gaussian distributed with variance a 2 , and the source uses 
the power P s P s .k in total to transmit the two symbols. The 
maximum source transmission rate at subcarrier k in the direct 
mode can be derived as 

Rk,2 = 2log 2 (l + P s>k ^^j bpts, (7) 

and it is achieved when the source uses the power s 2 5 ' to 
transmit each symbol. 

Note that we assume the same subcarrier is used by the 
source and the relays for transmitting a symbol in the relay- 
aided mode. In fact, optimized subcarrier pairing could also 
be implemented, which would further increase the degrees of 
freedom for optimization |[T6l . However, it would be more 
difficult to solve the RA problem, and that is why in the current 
work, the RA algorithms are designed under the assumption 
of nonoptimized subcarrier pairing. We will nevertheless use 
the insights gained here to guide future work, which will also 
consider subcarrier pairing. 

B. Formulation of the RA problem 

We consider the RA problem to maximize the sum rate, 
by optimizing the transmission mode and source power al- 
location for each subcarrier, as well as assisting relays and 
power allocation to them for every relay-aided subcarrier. To 
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formulate this RA problem, we define a binary variable tk 
which indicates the transmission at subcarrier k is in the relay- 
aided mode (resp. the direct mode) if tk = 1 (resp. tk = 0). 
Mathematically, the RA problem is formulated as 

K 

max ^ {tkRk,i + (1 - tk)Rk,2) 
fc=i 

K K 

s.t. ^P 8 , fc <l, 53P rtlfc <l > Vr i 6*, (8) 

k=l k=l 

P, fe >0,Vfc, Pr i ,fc>0,Vfc,Vr i G* J 
i fc e{0,l},Vfc, b k e{l,--- ,N},Vk, 

where {tk, bk, P s ,k, Prj,fc|Vrj € ^J/jVfc} are to be optimized 
by the RA algorithm. 

To facilitate analysis in the following sections, (0 can be 
formulated into the following equivalent form 



max /(x) 

X 

s.t. x G T> x , 
S(x) ^ 1, 



(9) 



where x represents the vector stacking all optimization vari- 
ables, T> x stands for the definition domain of x, /(x) denotes 
the sum rate, gr(x) = [51 (x), ■ ■ • ,g N+1 (x)} T where g x (x) = 
Efc=i p s,fc and i+ i(x) = £ fe=1 Pr,,k (* = !,••• > ^0 sta cks 
the sum power of the source and those of the relays, 1 
represents an (JV + 1) x 1 vector with every entry equal to 
1, and =^ denotes the entrywise "smaller than" inequality. 

III. The duality based RA algorithm 
The Lagrangian for ((9) is defined as 

L(x, / x) = /(x)+/x T (l-.g(x)), (10) 

where fi = [fi s ,fi Tl , • • • ,/j, In ] t , with fi a and fi Ti representing 
the dual variables related to the sum power constraints of the 
source and r,, respectively. The dual function is defined as 
d(fi) = max xg i) x i(x, ft), and we denote a x that maximizes 
L(x, ft) as x M . If there exist multiple x's that maximize 
L(x, fi), d(fi) is not differentiable at fi, and a subgradient of 
d(fi) at fi is equal to 1 — g(x M ), where any x that maximizes 
L(x, fi) can be chosen as x M l26l . 

It is important to note that the optimum objective value of 
(O, denoted by /*, always satisfies f*<d(fi) for any fi >p 0, 
where represents an (iV+1) x 1 vector with every entry equal 
to 0, and ^= denotes the entrywise "greater than" inequality 
11251 . The duality gap is defined as min M:/x ^o d(fi) — /*. Note 
that P x is not a convex set, and therefore the duality gap of 
(O might not be zero since the Slater constraint qualification 
is not satisfied l26l . Interestingly, it is shown in Appendix 
lAl based on the same idea proposed in ||231 that, when K is 
sufficiently large <j9j has zero duality gap. 

Assuming sufficiently large K, a duality based RA algo- 
rithm will be developed in this section to find a globally 
optimum solution to (O. This algorithm relies on finding the 
optimum dual variable fi* that minimizes d(fi) based on the 
subgradient method, i.e., updating fi with fi = [fi — S q (l — 
g(x fJ _))} + , where [fi] + represents a vector built from fi by 



only raising its negative entries to zero, and 5 q represents 
the step size used in the q-th iteration, until g(x„) =<! 1 and 
fi T (l — g(x M )) = are fulfilled. Then, x^* is a globally 
optimum solution to (0 as explained in Appendix A. If 5 q 
satisfies the diminishing conditions (i.e. lim g ^ +00 8 q — s- and 



2 0=1 ^1 = t 1 approaches fi* as q increases 

In practice, it may take an unaffordably large number 
of iterations before these optimality conditions are satisfied. 
To address this issue, the iteration can be terminated when 
g(x M ) =<; 1 and fi T (l — g(x^)) < e, where e is a prescribed 
small positive value. In this case, /(x M ) > /* — e follows 
because of 

/* - /(x M ) < d{fi) - /(x„) 

= A» T (1 - < e, 

which means that x M is a good approximation of a globally 
optimum solution to (|9j. 

The overall duality based RA algorithm is summarized 
in Algorithm [T] Specifically, fi is initialized with 1 at the 
beginning, and x M is found as the optimum solution to the 
Lagrangian maximization problem 



max L(x, ft) 



K 

E 

fe=i 



{t k Rkj + (l-t k )Rk,2) + 



K 



N 



K 



m s (i - E p -.*) + E ^ C 1 - E p ^)' 

k=l i=l fc=l 

s.t. P s , fe >0,Vfc, P r4)fc >0,Vfe,Vr < e*, (11) 
i fe e{0,l},Vfc, 6 fc G{l,..- > JV},Vfc > 



with Algorithm [2] developed in Section IIII-AI Furthermore, 
we choose S q = (Qi is a prescribed positive integer) 

which satisfies the aforementioned diminishing conditions. 

Algorithm 1 The duality based RA algorithm 

q = 1, fi = 1; {1 is a vector with every entry equal to 1} 
repeat 

q = q + U 

Find x^ with Algorithm |2j 

until g(x M ) =<! 1 and /x T (l — ff(x M )) < e 

x^ produced at last is output as an optimum solution to 



A. Algorithm to solve the Lagrangian maximization problem 

We can see that the optimum solution to ( TTTb 
can be found on a per subcarrier basis, i.e., 
{tfc(/i),6fc(M) ) -Ps,fc(M)>- p r 4 ,fc(M)|Vrj G *} representing 
the {tfe,6fc,P s ,fc,P ri ,fc|Vri G *} contained in x M for 
subcarrier k, can be found as an optimum solution to 



max 



Lk — tfePfe,! + (1 — tk)Rk,2 — MsPs.A 



N 



A^r; P"a,& 



i=l 

s.t. P s , fc >0, P ri>fc >0,V ri G (12) 

6fe G {1, ■ ■ • , iV}, t fe e{0,l}. 

When tk and 6^ are fixed, (fT2l is reduced to a con- 
vex optimization problem. When tk — 0, the maximum 
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L k is not influenced by b k , while this not the case when 
t k = 1 and bk is fixed. Based on the above analysis, an 
exhaustive-search based algorithm is designed to solve ([T2| >. 
To facilitate algorithm design, let's denote the maximum L k 
when t k = by Lk,o, the maximum Lk and the optimum 
{P B ,k, -P ri ,fc|V r 4 g *}whent fc = land b k is fixed by L kA (b k ) 
and {P s ,fc(&fe),P r ,,fe(&fe)|V r 4 g respectively. Note that 
L k ,i(b k ) and {P 8 ,fc(6 fe ), P r »,fc(&fe)|V r 4 € can be evaluated 
with Algorithm [3] developed in Section III.B. 

It is important to note that L k ,i(b k ) < L kl o if 
G s,v b .(k)(k) < G s ,d(fc), because of 



Lk,i(b k ) = log 2 (l +min(^,P s ^(6 fc )G s , rbfc(fc) (fc))) + X bk 
< log 2 (l + P s , fc (& fe )G s , d (fc)) + X bk (13) 
G s , d (fc)\ 



according to the KKT conditions that 



< 21og 2 1 + P s , k (h)- 



< L 



k,0, 



where X bk = -p s P s , k {b k ) - £) i=1 A* ri Prj,fc(&fc), and 
rjj. represents the % corresponding to ^fc(6fc) and 
{P,,fc(6fc),Pp i ,fc(6 fc )|V r,e *}. 



Algorithm 2 Algorithm to solve (fTTT i when /x is fixed 
for fc = 1 to K do 

if G s ,d(fe) > max ri g* G rii a(fc) then 
**(/*) = 0; 

P s>k (fj,) and P rii fc(/x),V rj 6 $ are evaluated by dl4t and 
(JT3J, respectively; 
else 

Compute Lfc.o with d 1 6b ; 

Find i k as the minimum i satisfying G sr . > G s ,d(fc); 

for 6 fe = i k to N do 

Find {P s ,k(&k)>Pr«>fc( & k)l v ^ G *} with Algorithm H 
when tfc = 1 and 6^ is fixed, then compute the L k> i(b k )', 

end for 

if L fe .o > max bkeBk L M (b fc ) then 
tfc(/i) = 0; 

P s ,fc(/x) and P Titk (fi), V rj £ ^ are evaluated by d!4t and 
| |15I >. respectively; 
else 
tfc(/*) = 1; 

b fc (/x) = argmax 6fceefc L k ,i(b k ); 

{P s ,fe(6fe),P ri ,fc(&fc)|V rj G *} when & fe = b k (fi) is 
assigned to {P Sjk (fi), P Ii:k (n)\V n G *}. 
end if 
end if 
end for 



Based on the above analysis, t k (fi), b k {fi), P Sik (fi), and 
{P ri) fc(/x)|V rj g can be found with one of the following 
procedures: 

1) when max rj G SyIi (k) < G Sj d(fc), t k (fi) = 0, because 
L k o > L k l (b k ) holds for every feasible value of 
6^ according to (fT~3~b . In this case, it can be derived 



log 2 (e) 



1 



P r( , fc (/i) = 0,Vr i g tf, 

"log 2 (e) 



L fc ,o = 21og 2 1 + G s , d (k) 



1 



Ms 



G Sjd (fc) 



(14) 
(15) 



(16) 



2) when max rj G SyYi (k) > G Sjd (fc), tfc(/i) can be deter- 
mined by an exhaustive-search based method, i.e., if 



max b fc 6{l,---,W} 



Lk,i(bk) > ife,o> <k(/f) = 1, otherwise 



tfe(/i) = 0. Let's denote i k as the minimum i satisfying 
G s , r ,(fc)(fc) > G sA (k). Note that when 1 < b k < i k -1, 
G s ,r bfc (fc)(fc) < G s ,d(fc), hence L k ,i(b k ) < L kfi holds 
for sure according to (1T3V This means that the compari- 
son of L k ,o with max ; , fce { 1 ... jv} L k} i(b k ) is equivalent 
to comparing L k ,o with max^^ L k>1 (b k ) where £ fc = 
{£&, • • • , TV}. Based on this idea, L k ,o is evaluated with 
( HoV and L kt x{b k ) is computed for all values of 6^ g 
6 fc with Algorithm |3] If L kfi > ma,x bkeBk L kA (b k ), 
t k (n) = 0, and P s , k (fi) and {P ri ,fe(/i)|V g are 
computed with (TT~4b and ( fl5l l. respectively. Otherwise, 
tfc(/i) = 1, &fc(/x) = argmax 6fcee(i L k ,i(b k ), and the 
{P s ,fc(6fc),P rj ,fc(6fe)|V rj g *} when 6 fc = b k (p) is 
taken as {P Sjfe (/ti), P ri , fc (/x)|V r t € 
In summary, the overall algorithm of finding x M is summa- 
rized in Algorithm[2] We will proceed with developing Algo- 
rithm Sto find L ktl (b k ) and {P s , k (b k ), P r „ fe (6 fc )|V n G 
in the next section. 

B. Algorithm to solve (1121 l w/ie« = 1 ant/ 6^ ;s ^xea" 

When t k = 1 and 6^ has a fixed value in (fT2l is 
equivalent to 

TV 

max L k = R kt i - /i s P s , fe - V" /i^Pr^fc 



/ / / 

i=l 



= log 2 (l + 7fc) - ^sPs.fc - / j Mr.Prj.fc 

s.t. 7 fc < Ps,fcG 8 , r .(*)(*), (17) 



N 

E- 

i=l 



7fe < P Sl fcG s , d (fc) + ( J2 \/ p ri,kG IiA (k)) , 



rje*fe(6 fc ) 

P rj , fe > 0. v r, g *, 

where 7^ is an intermediate optimization variable to guarantee 
the equivalence. 

It can readily be shown that ( fTTI i is a convex optimization 
problem. To solve it, one may formulate a set of equa- 
tions based on the KKT conditions and then solve them 
for {P s ,fc(6fc), P rii fc(6fc)|V ii G Vf}. It is very important to 
note that this method is effective only when the objective 
function and all the constraint functions are differentiable at 
the optimum solution. However, the second term in the right 
hand side of the second constraint in ( fTTI i is not differentiable 
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at P r ^ k = 0, Vr; G ffcOfc). This means that if 3 r, G * fe (6fe), 
Pii,k{bk) = 0, which might happen as shown later, the 
KKT conditions based method is not capable of finding that 
optimum solution. 

To address this issue, we solve ( fTTI i based on the idea that 
Ps,k(pk) is the optimum solution to 



max log 2 (l + 7fc) - n s P s ,k 
s.t. 7 fe < P s ,fcG Sjl . bfc(fc )(fc), 
Ik < Ps.kGsAk) + x, 

and {P ri .k(bk)\y r? G ^} is the optimum solution to 

N 



(18) 



max y^(-ju ri P ri , fe ) 



i=l 



S - t - E \J P n,kG Tit d(k) = y/x, 
f ri ,ft > 0,V ri G tf, 



(19) 



when x = x Q with 



51 J PnMGr^k) 

Specifically, x is first determined, and then (fT8l and (fT9l 
with x = x Q are solved to compute {P St k(b k ), Pr;,fc(&fc)|V r% G 
At first glance, this method seems confronted with a 
chicken-and-egg dilemma: though {P s , fc (& fe ), P u ,k{b k )\i r, G 
"f} can be computed by solving ( fT8l and ( fT9l once x Q is 
known, it seems that {P s , -Fri,fc(&fc)IV rj G \&} needs to 
be known first in order to compute x . In fact, this dilemma 
can be elegantly circumvented by using the sensitivity analysis 
in convex optimization theory to first determine x without 
knowing {P Syk (bk), Pr 4 ,fc(&fc)|V r$ G "J}, as elaborated in the 
following. 

1) Solutions to (fT8l and (fl9l l g/ven a;: Let's denote the 
optimum objective values of dT8l and (fl9l l by /i(x) and f2(x), 
respectively. Obviously, (fT8l is a convex optimization problem. 
Let's denote the optimum 7^ and the optimum dual variables 
associated with the first and second constraints of (fT8l by 
7fc(x), «fc(x), and (3 k (x), respectively. According to the KKT 
conditions of ( fT8l . 



Ms = G s . r ( fc) (fc)a fc (x) + C7 Sid (fc)/3fc(x) 



and 



log 2 e 



i + 



- 1 



(20) 



(21) 



_afc(x) + /3fe(x) 

should be satisfied. 

As for (fT9l , it can readily be derived that 
. when V r, G #k(&k), > 0, / 2 (x) = and the 

optimum P Ti , k to ( |19] l is 



if ri 

ifriG* fc (6 fc ). (22) 



Pruk = 

where 6fe = Er,e* fc (6 fc ) G r„d(*)//*r 4 



• when 3 r.; G ^ k {b k ), Mr, = 0. •M^) = 0, and the 
optimum P ri ,fc to (fT9l is equal to if ^ ^fc(frfc) or 
if rj G x I / fe(&fc) with [i Ti > 0. The optimum {P ri ,fc|V ij G 
^kipk),^ = 0} is any set of nonnegative values 
satisfying 



ri:riG*fc(6 h ),At r .=0 



Pn,kG n , d (k) 



(23) 



Moreover, it can readily be shown based on the Schwartz 
inequality that, V r; : r, G ^k(b k ) and = 0, the 
optimum P H . k satisfying d23l and minimizing the sum 

power of relays is P r k = jy '- G , k)r i - 

2) Finding x based on the sensitivity analysis: Let's de- 
note the Lk in dTTb computed with the optimum {7;-, P s , k } to 
(fT8l and the optimum {P ri ,fc|V r 4 G *} to £[9]l by L kil (x, b k ) 
when x > 0. Obviously, Lk,i(x,bk) = fi{x) + f2{x) and 
Lk,i(x,bk) < Lk,i(bk), since L kt i(b k ) is defined in Section 
III. A as the maximum for ( fTTI i, while L ki \(x,b k ) is 
the Lfc computed with the above mentioned {7fc,P s ,fc} and 
{P-i,fc|V rj G which are feasible for (fTTt . When x = x Q , 
the P s .fc and {P ri) fc|V r, G used for computing L k ,i(x, b k ) 
are equal to P s ,k(bk) and {P ri ,fc(&fe)|V r,; G *!/}, respectively, 
and therefore Lk.i(x ,b k ) = L k) i(b k ). This means that 
x = argmax^^o L k>1 (x, 6fe). 

To determine x Q , let's consider LI 1 (x, b k ) which represents 
the first order derivative of L k x(x,b k ) with respect to x. 
According to convex optimization theory, (3 k (x) represents the 
sensitivity of fi(x) with respect to x, i.e., f{(x) = /3 k (x) 
(please refer to pages 249-253 in j25l for more details). 
Therefore, 



L k,i(x,b k ) = 



P k (x) - l/9 bk 
(3 k (x) 



if Vr, G *fc(6 fc ),Ai r 
if 3t, G ^ k (b k ),fx r 



>0, 
= 0. 
(24) 



Based on the above analysis, the determination of x Q , 
/3fc(x ), afe(xo) and P s ,k(b k ) falls into one of the following 
cases: 

. when V n G * fc (6 fc ), /i r » > and G ^ (fc) < x D = 0. 
This is because L ky i(x, b k ) is a nonincreasing function of 
x since j3 k (x) G [0, ^ s /C7 Sjd (fc)] and L' k l (x,b k ) < for 
any x > 0. In this case, the first constraint in (TT8~t when 
relaxed, whereas the second one is saturated, 
since G S)d (fc) < Gs.r^ (*)(&)■ Therefore, a fc (x Q ) = 0, 



^^ ) = G^7(fe)' 7fe(^o) 



G 3 , d (fc)l 



1 



, and 



lk{x ) 
G s , d (fc) 



log 2 g 
Ms G s , d (fc) 



(25) 



> 



when V r 4 G *fc(6fc), /i ri > and g^f^y 
fik(x ) = since i' fe 1 (x Q , 6fe) = should be satisfied. 

M 5 -G s , d (fc)/ ebfc > Q because ^ 



,)(*) 



In this case, afc(x ) = 

is satisfied, and 7fe(x Q ) can be computed with (f2Tb . This 
means that both constraints in ( fT8l when 
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saturated, and therefore 



P s ,k{bk) 



Gs,r bk (k)(k) 

log 2 e 



1 



(26) 



and 



x a =P s , k (b k )A bk 
log 2 e 



fx s /A bk + l/e bk G s , r (fe) (fc) 



(27) 



where A bk = G s ^ bk{k) (k) - G B , d (k). 

. when 3 G ®k(h), Mr, = 0, /3 fc (x ) 

since Lj-. 1 (a; , b k ) = should be satisfied. In 

this case, a fc (:Co) = s — ^7— tct, and 7 fe (x ) = 
1 + 



r bfc (fc)(fc)log 2 e ^ 



This means that the first con- 



straint in (fl~8T > when X - — - X q IS saturated, whereas the 
second one is relaxed. Therefore, 



lki.X ) 

G s>1 (fc) (k) 



1 



log 2 e 

Ms G s>r (fc) (fc) 



1 - 



(28) 



and x- can be any value satisfying x > jk (x ) 
P s ,k( b k)G S:d {k) = x t h where 



1 



log 2 e 

Ms G s , r (fc) (fc) 



(29) 



Algorithm 3 Algorithm to solve (fTTT i when /x is fixed 



if V r, G *fe(6 fc ), Mr, > and G ^ (k) < -j- then 

P s fe(^fc) is evaluated by ( 1251 ); 
PrUfoO =0. Vr 4 G *; 
else if V n G 9 k (b h ), > and ^pry > ^- then 
Ps k(bk) is evaluated by d26t ; 

Vr.^f*), P H ,k(bk) = 0; 

V r; G P Ti ,k{bk) is evaluated by d22b with a; equal to 
a; computed with {27}; 

else if 3 ri G *fc(fefe), /U rj = then 
-P s k (bk) is evaluated by l !28t . 

Vr^* fc (6fc).i , r i ,fc(6fc) = 0; 

V n G *fe(6fc) with fj, ti > 0, P n ,k(b k ) = 0; 

V r t G *fe(6fe) with /Lt r< = 0, P H ,k(bk) = 

G r ,,d(*0*th 2 



end if 



ri ,d(fe)) 2 



After knowing x , we can find the optimum {P ri ,fe|V g 
to ([T9ll when x = x a as {F ri .fe(fofe)|V G Note 
that in the third case x can be any value no smaller than 
xth, and {-P ri ,fc(6fe)|V r< G *fc(6fc), /x ri = 0} can be any 
set of nonnegative values satisfying (l23l l with a; = ir . To 
improve the system energy efficiency, we choose x a = x t h, 
and V Ti : r< G * fc (6 fc ) and ^ = 0, P Iit k(h) = 
— ttxyt to minimize the sum power of 



(5Z ri :r ie 4. fc (6 fe ), Mr .=0 ' 

the relays. 



The overall algorithm of finding {P s> k(b k ), P Tt ,k( b k)\V T i G 
^17} is summarized in Algorithm [3] Based on the above 
analysis, it can be seen that the r\ k corresponding to 

{P s ,k( b k),Pn,k(bk)\V rj G and b k must be equal to or 
smaller than P s .k{bk)G s ^ bk{k) (k). 

IV. The iterative RA algorithm 

In case ^ has a nonzero duality gap, Algorithm Q] fails 
to find a globally optimum solution. To address this issue, 
we will develop in this section a coordinate-ascent based 
iterative algorithm which is suboptimum but always applicable 
regardless of the duality gap of (0. 

A. The iterative RA algorithm 

First of all, it should be noted that R k ,x < R k ,2 al- 
ways holds independently of {P s ,fcj -Pr;,fc|Vrj g \&} when 
b k satisfies G slb < G Sj d(fc). This means that when 

max rj G s>Ti (k) < G s .a(k), t k = is the optimum inde- 
pendently of b k and {P s ,k, -Pr;,fc|Vrj g ^/}, i.e., the direct 
transmission mode should always be used for subcarrier k. 

To simplify algorithm design, we assume t k is fixed as 
for every k g T) a = {fc| max rj G s>Ti (k) < G s>d (k)}. If 
k $l T) a , the optimum t k and b k must lie in the set {0, 1} x B k , 
where x is the Cartesian product operator. Furthermore, when 
{tk,b k \\/k $ T>a] is fixed, (0 is reduced to a concave 
maximization problem, which has zero duality gap since the 
Slater constraint qualification is satisfied. Thus a globally 
optimum {P Stk , P ri! fc|Vr. ; G \]/,Vfc} and the maximum sum 
rate when {t k ,b k \\/k ^ 2? a } is fixed can be found with a 
duality based algorithm as shown later. 

Now the difficulty of solving ^ lies in finding the optimum 
{tk,b k \Vk T> a }. To this end, one may use an exhaustive- 
search based algorithm. Specifically, this algorithm finds the 
maximum sum rate for each possible {t k ,b k \Vk ^ T> a } in 
Ilfe^x) ({A 1} x &k)- Then it chooses the best one as the 
optimum {t k ,b k \\/k ^ T) a }. However, the complexity of 
exhaustive search might be unaffordable for practical systems 
using a big number of subcarriers. 

To address this issue, we develop a coordinate-ascent based 
iterative algorithm which produces a successive set of RAs 
with nondecreasing sum rate until convergence. In the follow- 
ing, a superscript m added to a variable indicates that variable 
is produced at the m-th iteration to facilitate description. At 
the beginning, V k £ T) a , t\ and b\ are initialized as 1 
and N, respectively, i.e. every subcarrier not in T> a is set in 
the relay-aided mode with only the relay having the highest 
source-relay channel gain enabled for assisting. In the m-th 
iteration, {P s ™,P™ fc |Vr t G *,Vfc}, which is the optimum 
solution to ® when t k = tf and b k = b™, V k g V a , is 
first found with a duality based algorithm, namely Algorithm 
H] developed in Section IIV-BI This step can be interpreted 
as finding the optimum source/relay power allocation when 
{t k ,b k \Vk i V a } is fixed as {t?,V£\Vk i V a }. 

Then, {t k +1 ,b k +1 \M k (f. V a } which maximizes the 
sum rate is found when the power allocation is fixed 
as {P™ k ,P™ k \\/Ti G *,Vfc}. This step can be inter- 
preted as finding the optimum mode and assisting relays 
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Algorithm 4 The iterative RA algorithm 

V k € V a , t k = 0; V k £ V a , t\ = 1 and b\ = N; 

m = m + 1; 

repeat 

Find {-P s m fc , P™ fc |Vr.; G *,Vfc} as the optimum solution to l[8} 
when t fc = t^'and b k = b£\ V k £ V a with Algorithm |H 
for every k £ V a do 
if = then 



else 



m+l 



0; 



Compute C(fe"). flSliCb?). and i?£ 2 ; 
ifflr 1 (6H<-Rr 2 then 



= 0; 



else 



b? +1 



u k is set as the minimum i satisfying 77™ (6™) < 
-P s ™ C Siri (fc)(fc); 
end if 
end if 
end for 

771 = TO + 1; 

until Mk£V a , {t™ +1 , 6™ +1 } is equal to {t£\ 
ra.i^k.Vrj e *,Vfc} and {t^ft^.v* g v a ,t k ,b k y k e 
D a } produced in the last iteration are the suboptimum solution. 



when the source/relay power allocation is prescribed by 
{P s ™ , P r m fc |Vr, e Vfc}. Note that this can be accomplished 
on a per subcarrier basis, i.e., for every subcarrier k ^ T> a , 
t™ +1 and b™ +l are found to maximize the rate when the 
power allocation is fixed as {PJ^, P r TO fc|Vrj <G In this 

pm 

case, the rate is denoted by P™ 2 = 21og 2 (l + -^G Sy d(k)) 
if the direct mode is used. If the relay-aided mode with a 
given b k is used, the rate is denoted by R r kl (b k ) = log 2 (l + 
min(77^(6 fc ),P s m fe G sa . bfc(A . ) (fc))) where 77^ (6* ) represents the 
SNR of MRC at the destination, 



Let's first consider the evaluation of t 



m+l 



and when 



t™ = 0, meaning that subcarrier k was set in the direct mode 
when {P s ™ , -P r m fe |Vri e V&jVfc} was evaluated. Suppose the 
relay-aided mode with any b k l+1 is now used instead, the rate 
is reduced independently of because of 

KM +1 ) = lo g2(! + ^(Vk(bl l+1 ),P s n XG s ^ +lik) (k))) 

<log 2 (l + C(^" +1 )) 

= log 2 (l + P™ k G s4 (k))<R% 2 , 

where the equality in the third line is because of V r,; £ <]/, 
P™ fc = as will be shown in Section IIV-BI In order to 
maximize the rate, = 0, i.e. subcarrier k remains in 

the direct mode. 

Next, we consider the evaluation of t 7 k +1 and when 
t r k = 1, meaning that subcarrier k was set in the relay- 
aided mode with b k = bf when {P s m fc ,P " fc |Vr t e *,Vfc} 
was evaluated. As will be shown in Section lIV-Bl rjj^Q)™) < 
P™ k G s ^ bm {k){k) always holds and thus Rf tl (b%) = log 2 (l + 
Wkfik))- 11 can rea dily be seen from d5) that if b k > b™, 
Vi?(h) < vTK), while r,?(b k ) = Vf(bT) ^ h < bf, 
since V ^ ^.(fe™), P" L k — as will be shown in Section 
HV-Bl This means that V b k , P™i(&fc) < log 2 (l + rfkih)) < 
R k \{b k l ), i.e., if subcarrier k remains in the relay-aided mode, 
the rate can not be increased no matter which value is assigned 



to b k . In this case, t k n+1 and b™' +1 are determined in one of 
the following cases: 

• if R l n 2 > R k,-i( b k n )' C +1 is set as since the rate is 
increased if the direct mode is used; 

. if P™ 2 < P^O™), f k n+1 is still set as 1. b k n+1 is as- 
signed as the smallest b k satisfying R k n 1 (b k ) = P™i(&™)- 
Obviously, < b%, and ^ k (b k l+i ) is the biggest set 

of assisting relays that leads to the same rate as ^ k (¥ k ). 
The motivation behind this assignment is twofold. One is 
to guarantee the rate when b k 



is no smaller than 
when b k = frj™. The other is to increase the degrees of 
freedom for optimizing the power allocation in the next 
iteration, since V rj : r, ; g * fc (6™ +1 ) and r, ^ 
P r TO fc = always holds but P ruk is free to be optimized 
in the (to + l)-th iteration. It can easily be shown that 
is actually equal to the minimum i satisfying 
riTK) < P s m k G sMk) (k). Note that G SiI ^ +l(k) (k) > 

> G sA {k), which means that <= B k . 

It can easily be seen that a successive set of RAs with 
nondecreasing sum rate are produced as the iteration proceeds. 
After the algorithm converges, the RA produced by the last 
iteration is output as a suboptimum solution. This RA is at 
least a locally optimum solution to ^ (please refer to pages 
160-162 of ll26l for more details). In summary, the iterative 
algorithm is described in Algorithm [4] 

B. Algorithm to solve I© when t k and b k are fixed 

When tk = t™ and b k = &£*, V k V a , and t k = 0, 
V k e V a , © is equivalent to 

K 

max Y. Wk R *A + C 1 - ** J^a) + E R ^ 

k=l,k£V a keV a 
K K 

s.t. 5^P B , fc <l, ^P ri ,fe<l, Vr 4 , (30) 

k=l fe=l 

P s , fc >0,Vfc, Pr i ,fc>0,VA!,Vri6* I 

which has zero duality gap as shown earlier, and therefore 
can be solved with a duality based algorithm. To this end, we 
define the dual variables v s and v Ti related to the sum power 
constraints of the source and rj, respectively. Let's stack v B and 
-jVrJVr, e *} into a vector u, and {P s>k , P Iit k\ Vr, € *,Vfc} 
into a vector P. Then, we define the Lagrangian as 

K 

w(p, v)= y, (CPfe.i + (i - c)p fc , 2 ) + Yl Rk >* 

k=l,k$V a 



N 



N 



kev a 

K 



E v ^ - E + E ^ C 1 - E ■ 



k=l 



k=l 



The duality based algorithm looks for the optimum u with 
the subgradient based method. In each iteration, the optimum 
P that maximizes VF(P, u) subject to the constraints P s>k > 0, 
Pii,k > 0, V k, V ^ 6 denoted by P^, is found when 
v is fixed. We denote P s>k and P Ti , k contained in P„ as 
P s ,fc(i / ) and P Vi:k (v), respectively. Specifically, for subcarrier 
k with t r k = 1 and b k = b' k , the problem of finding 
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Algorithm 5 The duality based RA algorithm to solve ([H) 
when t k = tf and b k = bf, V k £ V a , and t k = 0, V k G V a 



source(square), relays(circle), destinations(diamond) 



j = l,v = 
repeat 

v s = \v. 



02 (I 



Ef =1 ^, fc H)] + >Vr lG *; 



i = i + 1; 

for fc = 1 to AT do 

If fc £ V a ov k <£V a with = then 

Ps,k(v) and Pr i ,fe(i / ) ! V ri are found with d!4t and dl51 
after replacing fj, s with i/ 8 ; 

else if k V a with i™ = 1 then 

P S) fc(i/) and P Tit k{v), V r, are found with Algo- 
rithm fj] after replacing {bk, fJ. B , fit- |V r, 6 'J?} with 
{&£>s,^|Vr, e *}; 

end if 
end for 

until P„ is feasible for d30t and 

A' A 

*(i - + e - E p **.fcH) < 6 

fc=l r;g* fc = l 

P 8j fc(i/), and P ri ,fc(i/), Vri 6 'J/ produced in the last iteration 
are the optimum solution to $8$ when tk = t^ and fofc = 
V fc £ V a , and i fc = 0, V n £ £> Q . 



{P S) /s(i/),P ri) fe(i/)|V n G *} is equivalent to (|T7j with 6 fc = 
b™, /i s — u s , and fx Ti = v Ti , V r.j G "f, and therefore can be 
solved with Algorithm [3] after replacing {bk, /U 8 , /i ri |V G "J} 
with {b™, v s , v Vi |V r, £ "J}. It is important to note that 
the r]k corresponding to {P S) fc(i/), P ri ,fc(t')|V G ^} and 
6^ = 6™ must be equal to or smaller than P S) k{v)G SjTbm (k), 
according to the analysis at the end of Section III.B. For 
subcarrier k with t r k n = or in T> a , the problem of finding 
{Pa,k{v), Pri,kiyW r « € *} is equivalent to (O with t k = 0, 
/U s = i/ s , and /i r; = i/ ri , V r; E "P, and therefore can be 
solved with ( [Pil l and < fT3T > after replacing /i s with v s . Note 
that Vr, € P ri ,fc(i/) = if k G £> a or t£* = 0, and 
Vr 4 £ Pr/.'fcM = if = 1. 

The overall algorithm is summarized as Algorithm [5] 
The step size and the termination condition are chosen in 
the same way as in Algorithm [T] Specifically, we choose 
5j = j£qJ (Q2 is a prescribed positive integer) which satisfies 
the diminishing conditions as proposed in [26]. Based on the 
analysis in Section [TTll the finally produced V v corresponds to 
a sum rate no smaller than z* — e (2* represents the maximum 
objective of (f30]l). 



V. Numerical experiments 

For illustration purposes, we consider an OFDM transmis- 
sion system aided by 6 DF relays shown in Figure [3] The 
source is located at the origin, the destination is located at the 
coordinate (0,-15), and rj, i = 1, • ■ ■ , 6, is located at the 
coordinates (-6,-7), (-4,-7), (-2,-7), (2,-7), (4,-7), 
and (6, —7), respectively. Note that all the aforementioned 
coordinate-related values have the unit of meter. The parame- 
ters are set as a 2 = -50 dBW, e = 0.1, Qi = Q 2 = 50, and 
the duration of one time slot is 1 millisecond. In numerical 
experiments, we found that the convergence of Algorithm 
Q] was always observed when K > 256. Thus we choose 
K = 256. 



<D —5 



-10 



-15 



-20 



-6 -4 -2 2 4 6 8 

x-axis (meter) 

Fig. 3. A DF relays aided OFDM transmission system considered in 
numerical experiments. 



channel gain from the soure to a relay or the destination 




100 150 
subcarrier k 



channel gain from a relay to the destination 




100 150 
subcarrier k 



Fig. 4. A random channel gain realization. 



We generate every channel according to the following 
assumptions. First, it is modeled as a 6-tap delay line, and the 
average received power at a distance of d is equal to £,d~ 3 , 
where £ represents a log-normal shadowing effect of 1 dB 
(i.e., 101og 10 (£) is Gaussian distributed with zero mean and 
variance equal to 1 dB). Second, we assume the amplitude 
of the ?'-th tap is a circularly symmetric complex Gaussian 
random variable with zero mean and variance as of = <7oe _3 \ 
i = 0, 1, • • • ,5. 

In order to illustrate the effectiveness of the proposed algo- 
rithms, we consider a heuristic RA algorithm against which 
the proposed RA algorithms are compared. This heuristic 
algorithm allocates the source sum power uniformly to all 
subcarriers, i.e., P s<k — -h. Then, for every r,, O rj = 
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Fig. 5. For the RA evaluated by each algorithm, the sum rate and the TETIB 
when P s = Pr , , V r; e and P s varies from 10 to 50 dBW. 
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Fig. 6. The RA evaluated by Algorithm [T]for the random channel realization 
when P s = P n = 20 dBW, V r, 6 *. 



{k\k ^ T> a) G StTi (k) = max r .g$ G s ,r (k)} which contains 
every subcarrier not belonging to T> a and at which r< has the 
maximum source to relay channel gain, is first found. The 
sum power of r.; is uniformly allocated to all subcarriers in 
£l Ti (i.e., if k ^ Vt ti , P Iit k = 0, otherwise P n ,k = Td7T wnere 
|O ri | denotes the number of subcarriers in For every 

subcarrier k S T> a , the direct transmission mode is used. For 
every subcarrier k ^ V a , Rk } 2 is first computed, and then 
the Rk,\ when a single with the maximum G SyIi (k) assists 
relaying is computed. If Rk,2 is equal to or greater than the 
computed Rk,i, the direct mode is used. Otherwise, the relay 
mode is used, and only the used for computing Rk,x assists 
relaying. 

We have generated a single random channel realization as 
shown in Figure |U and tested the three RA algorithms when 
P s = P Ti , V n e and P s varies from 10 to 50 dBW. For 
the RA evaluated by each algorithm, the sum rate and the 
total energy of the source and the relays for transmitting per 
information bit (TETIB) are shown in Figure |5] We can see 
that for each P s , either Algorithm [T] or Algorithm Uproduces 
a better RA with a higher sum rate and a lower TETIB than 
the heuristic algorithm, and Algorithm [Uproduces a better RA 
than Algorithm S] Moreover, the RA evaluated by Algorithm 
1 or Algorithm 4 results in a much smaller TETIB than that 
computed by the heuristic algorithm especially when P s is 
low. For the sake of space limitation, the RAs evaluated by 
the three algorithms for the random channel realization only 
when P s = P u = 20 dBW, Vr t e are shown in Figure |6]|7] 
and|H] respectively. When the RA evaluated by Algorithm 1 is 
used, r3 and are enabled to assist relaying simultaneously at 
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a few subcarriers, while only 1-4 does when the RA computed 
by Algorithm 4 or the heuristic algorithm is used. 

We have also generated 100 random channel realizations, 
and tested the three RA algorithms when P s = P Xi , V £ 
and P s varies from 10 to 50 dBW. For the RAs evaluated by 
each algorithm, the average sum rate and the average TETIB 
are shown in Figure [9] It can be seen that for each P s , either 
Algorithm |4] or Algorithm Q] produces better RAs with a 
higher average sum rate and a lower average TETIB than the 
heuristic algorithm, and Algorithm[T]produces better RAs than 
Algorithm U Moreover, the RAs evaluated by Algorithm 1 
or Algorithm 4 result in a much smaller average TETIB than 
those computed by the heuristic algorithm especially when P s 
is low. 

VI. Conclusion 

We have considered the sum rate maximized RA problem 
in a point to point OFDM transmission system assisted by 
multiple DF relays subject to the individual sum power con- 
straints of the source and relays. In particular, one or several 
relays may cooperate with the source to transmit at every relay- 
aided subcarrier. We have proposed two RA algorithms which 
optimize the assignment of transmission mode and source 
power for every subcarrier, as well as the optimum assisting 
relays and the power allocation to them for every relay aided 
subcarrier. The effectiveness of the two algorithms has been 
illustrated by numerical experiments. 

Appendix A 

When the duality gap is equal to zero, the duality based 
algorithm is one of the state-of-the-art methods to find a 



globally optimum solution to (0. There exist two important 
theorems justifying the ability of the duality based algorithm to 
find a globally optimum solution to (0. The first (Proposition 
3.3.4 in 11261 ) is that x M is a globally optimum solution to 
©, if it is feasible to © and - .g(x M )) T = 0. The 
second (Proposition 5.1.4 in J26)) is that the satisfying 
the aforementioned conditions, denoted by fi*, must minimize 
d{n) when fi )p 0. 

In general, the duality gap of (0 can be studied with a 
visualization based method proposed in 11261 . Specifically, we 
can plot a cloud of points collected in the set S = {[p, w]\p = 
g(x),w = /(x),x £ 2? x } in the hyperplane of [p,w] shown 
in Figure [10] Most interestingly, d(/x) is equal to the in- 
coordinate of the highest intersection between the line p = 1 
and a line passing through S and perpendicular to the vector 
[—//., 1]. As shown in Fig. [10] the duality gap is equal to zero 
if the w-coordinate of the upper border of S is a concave 
function of p. Mathematically, a point on the upper border of 
S has the coordinate (p, /(x p )), where 

x p = arg max /(x). (31) 

x:xED x , p(x)— p 

Based on the same idea first proposed in ||231 . we can show 
/(x p ) is a concave function of p when K is very large. This 
is equivalent to show that for any 9 £ [0, 1], 

/(x epi+(1 _ e)p2 ) > 0/(x Pl ) + (1 - 0)/(x P2 ) (32) 

holds for any pi and P2 stacking the sum power of the source 
and relays. Note that the above condition can be interpreted 
in a very interesting way as follows. Let's adopt the RAs x pi 
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Fig. 10. The visualization of the duality gap in the hyperplane of [p, tu]. 

and x P2 in the 9 and 1 — 9 portions of the whole transmission 
duration, respectively, which is called 9 time sharing of x pi 
and x P2 hereafter. In this way, an average sum rate 6*/(x pl ) + 
(1 — f?)/(x P2 ) can be achieved with an average sum power 
(9pi + (1 — 0)p2- This means that showing the validity of (13 2 b 
for any 9 <= [0, 1], is equivalent to show that the optimum RA 
for the sum power #pi + (1 — (9)p2 provides a higher sum rate 
than 9 time sharing of x pi and x P2 . 

To this end, we show in the following that when K is 
very large, a RA x' € D x , which is of the sum power 
#pi + (1 — 9)y>2 and yields a sum rate equal to 9 time 
sharing of x pi and x P2 , can be constructed by 9 spectrum 
sharing of x pi and x P2 , i.e., taking the entries of x pi at 9 
portion of all subcarriers, and the entries of x P2 at all the 
other subcarriers to construct x'. More specifically, we first 
divide all subcarriers into S subbands, each consisting of N s 
adjacent subcarriers experiencing the same channel conditions. 
Suppose K is sufficiently large so that N s is also very large, 
and V 9 6 [0, 1] there exists an integer N$ in {1, • • • ,N S } 
satisfying 9 rj Since the subcarriers in each subband 
experience the same channel conditions, the entries of x pi 
(or x P2 ) at the subcarriers in the same subband should be the 
same. Let's construct x', which in every subband adopts for 
the first Ng subcarriers the same entries as in x pi , and for the 
remaining N s — Ng subcarriers the same entries as in x P2 . We 
can easily show that the conditions x' € P x and 

«?(x')«#Pl + (l-#) P 2, 

/M«9/tPi) + (l-«)/W 

are all satisfied. Thus /(x epi+(1 _ 0)p2 )>/(x') « 0/(pi) + 
(1 — #)/(pa) holds. Therefore, the duality gap of ((9) is equal 



to zero when K is very large. 
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